
FIELD OF THE INVENTION 

The object of the present invention is a method intended for gradual deformation of 
representations or reahzations, generated by sequential simulation, of a not necessarily 
Gaussian stochastic model of a heterogeneous medium, based on a gradual deformation 
algorithm of Gaussian stochastic models. 

The method according to the invention finds applications in underground zones 
modelling intended to generate representations showing how a certain physical quantity 
is distributed in an underground zone (permeability z for example) and best compatible 
with observed or measured data : geologic data, seismic records, measurements 
obtained in wells, notably measurements of the variation with time of the pressure and 
of the flow rate of fluids from a reservoir, etc. 

BACKGROUND OF THE INVENTION 

In patent application FR-98/09,018 is described a method intended for gradual 
deformation of a stochastic (Gaussian type or similar) model of a heterogeneous 
medium such as an underground zone, constrained by a set of parameters relative to the 
structure of the medium. This method comprises drawing a number p (p=2 for example) 
of realizations (or representations) independent of the model or of at least part of the 
selected model of the medium from all the possible realizations and one or more 
iterative stages of gradual deformation of the model by performing one or more 
successive linear combinations of p independent initial realizations, then composite 
realizations successively obtained possibly with new draws, etc., the coefficients of this 
combination being such that the sum of their squares is 1. 
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Gaussian or similar models are well-suited for modelling continuous quantity fields 
and they are therefore ill-suited for modelling zones crossed by fracture networks or 
channel systems for example. 

The most commonly used geostatistical simulation algorithms are those referred to 
5 as sequential simulation algorithms. Although they are particularly well-suited for 
simulation of Gaussian models, they do not imply in principle a limitation to this type of 
model. 

A geostatistical representation of an underground zone is formed for example by 
subdividing it by a network with N meshes and by determining a random vector with N 

10 dimensions Z = (Zi, Z2,...Zn) best corresponding to measurements or observations 
obtained on the zone. As shown for example by Johnson, M.E. ; in « Multivariate 
Statistical Simulation » ; Wiley & Sons, New York, 1987, this approach reduces the 
problem of the creation of an N-dimensional vector to a series of N one-dimensional 
problems. Such a random vector is neither necessarily multi-Gaussian nor stationary. 

15 Sequential simulation of Z first involves the definition of an order according to which 
the N elements (Zi, Z2,...Zn) of vector Z are generated one after the other. Apart from 
any particular case, it is assumed that the N elements of Z are generated in sequence 
from Zi to Zn. To draw a value of each element Zi, (i = 1, N), the following 
operations have to be carried out : 

20 a) building the distribution of Zi conditioned by (Zi, Z2„Z,.\) 



Fe (Zi) = P(Zi < Zi / Zi, Z2,. . .Zi.,) ; and 



b) drawing a value of Zi from distribution Fc (Zi). 




In geostatistical practice, sequential simulation is frequently used to generate multi- 
Gaussian vectors and non-Gaussian indicator vectors. The main function of sequential 
simulation is to determine conditional distributions Fc (Zi) (i = 1,..., N). Algorithms and 
softwares for estimating these distributions are for example described in : 

5 - Deutsch, C.V. et al, « GSLIB (Geostatistical Software Library) and User's Guide » ; 
Oxford University Press, New York, Oxford 1992. 

Concerning drawing the values from distribution Fc (zj), there also is a wide set of 
known algorithms. 

We consider the inverse distribution method by means of which a realization of Zi : 
10 Zi = Fc'^ (Ui) is obtained, where Ui is taken from a uniform distribution between 0 and 1. 
A realization of vector Z therefore corresponds to a realization of vector U whose 
elements Ui, U2, Un, are mutually independent and evenly distributed between 0 and 
1. 

It can be seen that a sequential simulation is an operation S which converts a 
15 uniform vector U = (Ui, U2, . . Un) to a structured vector Z = (Zi, Z2, . . ., Zn) : 

Z = S(U) (1). 

The problem of the constraint of a vector Z to various types of data can be solved 
by constraining conditional distributions Fc (Zi) (i = 1, N) and/or uniform vector U = 
(U„U2, ...,Un). 

20 Recent work on the sequential algorithm was focused on improving the estimation 

of conditional distributions F^ (Zi) by geologic data and seismic data integration. 
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An article by Zhu, H. et al can be mentioned for example : « Formatting and 
Integrating Soft Data : Stochastic Imaging via the Markov-Bayes Algorithm » in Soares, 
A., Ed. Geostatistics Troia 92, vol.1 : Kluvi^er Acad. Publ, Dordrecht, The Netherlands, 
pp.1-12, 1993. 



pressures from well tests and production records, unless a severe linearization is 
imposed. Furthermore, since any combination of imiform vectors U does not give a 
uniform vector, the method for gradual deformation of a stochastic model deve lopped in 
the aforementioned patent cannot be directly applied v/ithin the scope of the sequential 
1 0 technique reminded above . 

The method according to the invention thus allow^s to make the two approaches 
compatible, i.e. to extend the formalism developped in the aforementioned patent to 
gradual deformation of realizations, generated by sequential simulation, of a not 
necessarily Gaussian stochastic model. 



The method allows gradual deformation of a representation or realization, generated 
by sequential simulation, of a not necessarily Gaussian stochastic model of a physical 
quantity z in a heterogeneous medium such as an imdergroimd zone, in order to 
constrain it to a set of data collected in the medium by previous measurements and 
20 observations relative to the state or the structure thereof. 

It is characterized in that it comprises applying an algorithm of gradual deformation 
of a stochastic model to a Gaussian vector (Y) having a number N of mutually 
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However, this approach cannot be extended to integration of non-linear data such as 
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DEFINITION OF THE METHOD 
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independent variables that is connected to a uniform vector (U) with N mutxially 
independent uniform variables by a Gaussian distribution function (G), so as to define a 
chain of realizations u(t) of vector (U), and using these realizations u(t) to generate 
realizations z(t) of this physical quantity that are adjusted in relation to the (non- linear) 



According to a first embodiment, the chain of realizations u(t) of uniform vector (U) 
is defined from a linear combination of realizations of Gaussian vector (Y) v^hose 
combination coefficients are such that the sum of their squares is one. 

According to another embodiment, gradual deformation of a number n of parts of 
10 the model representative of the heterogeneous model is performed while preserving the 
continuity between these n parts of the model by subdividing uniform vector (U) into a 
number n of mutually independent subvectors. 

BRIEF DESCRIPTION OF THE DRAWINGS 

Other features and advantages of the method according to the invention will be 
15 clear from reading the description hereafter of a non limitative example, with reference 
to the accompanying drawings wherein : 

Figure 1 shows the medial layer of a realization of a facies model used as a 
reference, generated by sequential simulation of indicatrices, 

- Figure 2 shows the variation with time of the pressure obtfiined in a well test for the 
20 reference model. 

Figures 3 A to 3E respectively show five initial realizations of the medial layer of a 
reservoir zone, constrained only by the facies along the weU, 
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data. 



- Figures 4A to 4E respectively show, for these five realizations, the bottomhole 
pressure curves in the reference model compared with those corresponding to the initial 
models, 

- Figures 5 A to 5E respectively show five realizations of the medial layer of the facies 
model conditioned to the facies along the well and adjusted in relation to the pressure 
curve obtained by well tests, 

- Figures 6A to 6E respectively show, for the five realizations, the bottomhole 
pressure curves in the reference model compared with those corresponding to the 
adjusted models. 

Figures 7A to 7E respectively show how the objective fiinctions respectively 
corresponding to these five examples vary with the number of iterations, 

- Figures 8A to 8E show the gradual deformations generated by an anisotropy 
coefficient change on a three-facies model generated by sequential simulation of 
indicatrices, and 

Figures 9A to 9E show the local gradual deformations of a three-facies model, 
generated by sequential simulation of indicatrices. 

DETAILED DESCRIPTION OF THE METHOD 

We consider a study zone that is subdivided by an N-mesh grid and we try to build 
realizations or representations of a stochastic model of a certain physical quantity z 
representing for example the permeability of the formations in the zone. The wanted 
model must adjust to data obtained by measxirements or observations at a certain 
number of points, and notably to non- linear data. 
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Adjustment of a stochastic model to non-linear data can be expressed as an 
optimization problem, f**^ = (f.°^ f, f,-"^' ... f^'^'^^) designates the vector of the non- 
linear data observed or measured in the medium studied (the reservoir zone), and f = (fi, 
fs, fa . , . fp) the corresponding vector of the responses of the stochastic model of the 
permeability Z = (Zi, Z2,..., Zn). The problem of constraining stochastic model Z by 
observations consists in generating a realization z of Z which reduces to a rather low 
value an objective function that is defined as the sxmi of the weighted rms errors of the 
responses of the model in relation to the observations or measurements in the reservoir 
zone, i.e. : 

where coi represents the weight assigned to response 5. Functions (i=l, 2, p) and 
objective fiinction O are functions of vector Z. We are thus faced with an optimization 
problem of dimension N. 

In order to extend the formalism developped in the aforementioned patent to the 
gradual deformation of realizations generated by not necessarily Gaussian sequential 
simulation, we start from a Gaussian vector with N variables Yj with i = 1, 2, N, 
mutually independent, of zero mean and of variance equal to 1, and N mutually 
independent uniform variables Ui, U2, U3,...Un are defined by : 

Ui = G(YO Vi=l,2, ...,N 

where G represents the standardized Gaussian distribution function. 

Assuming this to be the case, the gradual deformation algorithm developped within 
a Gaussian fi-ame is applied to the Gaussian vector Y = (Y,, Y2,..., Yn) in order to build 





a continuous chain of realizations of uniform vector U = (Ui, U2,..., Un). Given two 
independent realizations Ya and yt of Y, the chain of realizations u(t) of vector U 
obtained with the following relation is defined : 

u(t) = G(y3 cos t + yb sin t) (2). 

5 For each t, u(t) is a realization of vector U. A vector z(t) which is, for each t, a 

realization of random vector Z is then obtained by sampling of the conditional 
distribution Fc (Zi) (i=l, 2,..., N) using the elements of vector u(t). Parameter t can 
consequently be adjusted as in the Gaussian case so as to adjust z(t) to non-linear data. 
This procedure is iterated until satisfactory adjustment is obtained. 

10 Adjustment of a facies model to pressure data obtained by means of well tests 

In order to illustrate application of the stochastic optimization method defined 
above, we try to adjust a stochastic reservoir model to pressure data obtained by means 
of well tests. Building of the reservoir model derives fi-om a real oil formation 
comprising three types of facies : two reservoir facies of good quaUty (facies 1 and 2) 
15 and a reservoir facies of very bad quality (facies 3). Table 1 defines the petrophysical 
properties of the three facies : 
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Facies 1 


10 


10 


10 


17 


2.1857 


Facies 2 


1 


1 


1 


14 


2.0003 


Facies 3 


0.1 


0.1 


0.001 


9 


1.8148 



In order to represent the specific facies distribution of the oil formation, a binary 
realization is first generated to represent facies 3 and its complement. Then, in the 
complementary part of facies 3, another binary realization independent of the first one is 




generated to represent facies 1 and 2. The formation is discretized by means of a regular 
grid pattern of 60x59x15 blocks 15mxl5mxl.5m in size. An exponential variogram 
model is used to estimate the conditional distributions. The main anisotropy direction is 
diagonal in relation to the grid pattem. The ranges of the variogram of facies 3 in the 
5 three anisotropy directions are 300m, 80m and 3m respectively. The ranges of the 
variogram of facies 1 and 2 in the three anisotropy directions are 150m, 40m and 1.5m 
respectively. The proportions of facies 1, 2, 3 are 6%, 16% and 78% respectively. 

A well test has been carried out by means of a finite-difference well test simulator 
as described by : 

10 Blanc, G. et al : « Building Geostatistical Models Constrained by Dynamic Data - A 
Posteriori Constraints » in SPE 35478, Proc. NPF/SPE European 3D Reservoir 
Modelling Conference, Stavenger, Norway, 1996. 

The medial layer of a realization used as the reference model for this validation can 
be seen in Figure 1. The section of the well that has been drilled runs horizontally 
15 through the medial layer of the reservoir model along axis x. The diameter of the well is 
7.85cm, the capacity of the well is zero and the skin factors of facies 1, 2 and 3 are 0, 3 
and 50 respectively. The synthetic well test lasts for 240 days with a constant flow rate 
of 5 m^/day so as to investigate nearly the entire oil field. Figure 2 shows the pressure 
variation with time. 

20 The objective was to build realizations of the facies model constrained by the facies 

encountered along the well and by the pressure curve obtained during well testing. The 
objective function is defined as the sum of the rms differences between the pressure 
responses of the reference model and the pressure responses of the realization. Since the 
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dynamic behaviour of the reservoir model is mainly controlled by the contrast between 
the reservoir facies of good and bad quality, the binary realization used to generate 
facies 1 and 2 has been fixed first and only the binary realization used to generate facies 
3 has been deformed for pressure data adjustment, 

5 The pressure responses resulting from the well tests for the five realizations of 

Figs.3 A to 3E are different from that of the reference model, as shown in Figs.4A to 4E. 
Starting respectively from these 5 independent realizations, by using the iterative 
adjustment method above, we obtain, after several iterations, five adjusted realizations 
(Figs.SA to 5E) for which the corresponding pressure curves are totally in accordance 

10 with those of the reference model, as shown in Figs.6A to 6E. 

Gradual deformation in relation to the structural parameters 

In many cases, sufficient data for deducing the structural parameters of the 
stochastic model : mean, variance, covariance function, etc, is not available. These 
structural parameters are often given in terms of a priori intervals or distributions. If 

15 their values are wrong, it is useless to seek a realization adjusted to non-linear data. It is 
therefore essential for applications to be able to perform a gradual deformation of a 
realization with simultaneous modification of random numbers and structural 
parameters. The sequential simulation algorithm defmed by equation (1) makes it 
possible to change, simultaneously or separately, structural operator S and uniform 

20 vector U. Figs. 8A to 8E show the gradual deformations obtained for a fixed realization 
of imiform vector U when the anisotropy coefficient is changed. 
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Local or regionalized gradual deformation 

When the observations are spread out over dififerent zones of a formation studied, 
an adjustment using global deformation would be ineffective because the accordance 
improvement obtained in a zone could deteriorate it in another zone. It is therefore 
preferable to apply gradual deformations zone by zone. Consider a subdivision of vector 
U into a certain number n of mutually independent subvectors IP,..., U°, which 
allows to perform their gradual deformation individually. Separate application of the 
gradual deformation algorithm to each subvector U*, IP,..., U" allows to obtain a 
function of dimension n of uniform vector U : 

G(rjcosr,+7^*sin/,) " 
G(r/cos/2+J;'sin t^) 

G(7; cos + sin o_ 

where Y'a and Y^b for any i = 1, 2,..., n, are independent Gaussian subvectors. For a set 
of realizations of and Y*b, a problem of optimization of n parameters ti, t2,..., t^ is 
solved to obtain a realization that improves or maintains the data adjustment. This 
procedure can be iterated until satisfactory adjustment is obtained. 

Gradual local deformations thus allow to significantly improve the adjustment 
speed in all the cases where measurements or observations are spread out over different 
zones of the medium. 

The effect of this gradual local deformation on the three-facies model of Figs. 9 A to 
9E can be clearly seen in these figures where only the delimited left lower part is 
affected. 
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The method accordmg to the invention can be readily generalized to gradual 
deformation of a representation or realization of any stochastic model since generation 
of a realization of such a stochastic model always comes down to generation of uniform 
numbers. 



